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Abstract 

We consider a model of the Universe in which the matter content is in the form of 
discrete islands, rather than a continuous fluid. In the appropriate limits the resulting 
large-scale dynamics approach those of a Friedmann-Robertson- Walker (FRW) universe. 
The optical properties of such a space-time, however, do not. This illustrates the fact that 
the optical and 'average' dynamical properties of a relativistic universe are not equivalent, 
and do not specify each other uniquely. We find the angular diameter distance, luminosity 
distance and redshifts that would be measured by observers in these space-times, using 
both analytic approximations and numerical simulations. While different from their coun- 
terparts in FRW, the effects found do not look like promising candidates to explain the 
observations usually attributed to the existence of Dark Energy. This incongruity with 
standard FRW cosmology is not due to the existence of any unexpectedly large structures 
or voids in the Universe, but only to the fact that the matter content of the Universe is 
not a continuous fluid. 
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1 Introduction 

One of the most basic assumptions made in modern cosmology is the idea that the geometry 
of space is homogeneous and isotropic, resulting in the line-element 



where k is the conformal curvature of space, and a(t) is the time-dependent scale factor of 
the Universe. The simplification that results from this assumption is remarkable: Einstein's 
equations, a set of 10 coupled, non-linear partial differential equations in 4 variables, reduce to 
a single ordinary differential equation in one variable: 



The over-dot here is a derivative with respect to time, and p is the (spatially constant) energy 
density of the continuous fluid that is assumed to permeate the whole of space. Such models 
are known as Friedmann- Robert son- Walker (FRW) cosmologies, and are ubiquitous. 

The unreasonable simplicity offered by the assumptions outlined above has allowed enormous 
progress to be made in understanding the cosmological models that result from them, and the 
ways in which they can be compared to astronomical observations. So great is this success that 
there now exists a growing sense that such models must be the only ones that are suitable to 
describe our Universe. These sentiments are bolstered by observations of the Cosmic Microwave 
Background (CMB), which is isotropic to within one part in a hundred thousand. 

However, albeit effective, the homogeneity and isotropy of space is still only an assumption, 
and the isotropy of the CMB, while consistent with a homogeneous and isotropic space-time, 
does not necessitate it P0|. Furthermore, it now seems to be the case that these assumptions 
lead inevitably, and unenviably, to the conclusion that the Universe should be filled with an 
exotic fluid that behaves repulsively under gravity. This fluid, dubbed Dark Energy, has to 
make up ~ 2/3 of the total energy budget of the Universe and should be responsible for driving 
a late period of apparently accelerating expansion [2j [3j HJ [5] . Unfortunately, the existence of 
such a fluid provides colossal theoretical challenges, and despite ongoing efforts to find it, it has 
yet to be detected directly. It therefore seems natural to retrace our steps, and to scrutinise 
more critically the assumptions that led us to deduce its existence in the first place. 

Inhomogeneous cosmology has been an active field of research for many decades now [5]. 
One of the simplest examples, which places us directly at the centre of the Universe, is the 
Lemaitre-Tolman-Bondi (LTB) model: A spherically symmetric, and dust dominated model. 
Once again it is implicitly assumed that there exists a continuous fluid that permeates all of 
space-time. It is also usual to assume that the spherical inhomogeneity in this model takes 
the form of a giant under-density at the centre of the Universe - i.e that we live in a void. 
LTB models have been extensively studied over the past few years, and their behaviour and 
optical properties are now becoming well understood [7]- [26]. They can mimic an FRW universe 
with Dark Energy, but are required to have a fairly intricate structure in order to be able to 
simultaneously match supernova observations and measurements of primary anisotropies in the 
CMB. 
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A natural next step is to consider multiple LTB or Schwarzschild space-time patches, em- 
bedded in an FRW background. The inhomogeneities that result are then part of a continuous, 
and otherwise constant fluid that constitutes an exact solution of Einstein's equations. Such a 
Universe could conceivably be treated as being homogeneous and isotropic on the largest scales, 
but is endowed with substantial inhomogeneity on smaller scales. Dubbed as "Swiss Cheese" , 
these models have been explored in great detail; from the ground breaking work of Kantowski 
[27J, up to the present day [2E]-|22]- Again, substantial effort has gone into understanding the 
optical properties of Swiss Cheese universes, and they have been shown to be broadly similar, 
but not identical, to those of FRW. 

Beyond the realm that can be easily investigated with exact solutions, one will also be in- 
terested in approximations that may allow further insight to be had. One scheme for including 
inhomogeneity without the aid of exact solutions was provided by Dyer and Roeder [331 G3]- 
They ignored the influence of shear on the evolution of bundles of null geodesies, and taking 
FRW values for the Hubble rate and redshift found a neat method of approximating observ- 
ables in an inhomogeneous universe. This method was generalised further by Mattson [35], who 
included the effect of an inhomogeneous Hubble rate. It was claimed in [35] that with an appro- 
priately changing H(z) the effect of inhomogeneities could entirely account for the deviations 
from Einstein-de Sitter (EdS) cosmology that are usually attributed to Dark Energy. However, 
a lack of any model on which to base these results, and the ad hoc way in which inhomogeneous 
expansion is dealt with, makes these claims appear somewhat speculative. 

A further avenue of research that has attracted considerable interest is the idea that the 
structures that form in the Universe could have a 'back-reaction' effect on the cosmological 
expansion [36]-|52j. The idea here is based on the fact that the operation of averaging the 
geometry of space does not, in general, commute with the operation of evolving a space-like 
surface forward in time. Such complications mean that when evolving forward an 'averaged' 
homogeneous and isotropic space one should include corrections to Einstein's equations (which 
are valid for the true geometry, and not the averaged one). Some have argued that the results 
of these corrections may be large enough to entirely account for Dark Energy, while others 
maintain the opposite position. These studies are considerably complicated by the implicit 
difficulties involved with averaging the geometry of space in General Relativity. 

In this paper we choose to take a different path, and completely break with FRW space-time. 
We want to work with a model of the Universe in which all of the matter content is in the form 
of discrete islands of mass, scattered in otherwise empty space. To a first approximation, this 
appears to be what we see on the night sky: Scattered points of light, organised into a variety 
of structures and meta-structures. Our Universe, at its most basic, is made up of galaxies with 
a typical mass roughly that of the Milky Way, dominated by dark matter and with a number 
density of ~ 1 galaxy per cubic mega-parsec. 

Such a universe should have distinctive properties. Light will propagate in the empty spaces 
between the islands, and will no longer pass through a continuous fluid. The geometry of these 
two situations is very different, and should be expected to result in different optical properties, 
even if the large-scale dynamics are equivalent [53]. To see why this is so, consider the form of 
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the Sachs optical equations [M] (which we will discuss in detail later on): 



— + 6 + o- 



R 




C. 



Here 6 and a are the expansion and shear scalars, respectively, C represents the Weyl tensor, 
and R the Ricci tensor. The angular diameter and luminosity distances are then given by 
integrals of 9. In a universe where light travels primarily through empty space, the driving 
terms in the equations above will be R = and The corresponding terms in a spatially 

flat FRW universe will be C = 0, due to its conformal equivalence to Minkowski space, and 
R 7^ 0, due to its continuous mass distribution. These two fundamentally different types of 
curvature can be seen to have correspondingly different effects on bundles of null geodesies, and 
hence these effects should be taken into account in observational cosmology. 

It is, or course, a very difficult proposition to calculate the geometry of space-time associated 
with arbitrarily scattered islands of matter. Fortunately, if we return the underlying principle 
behind FRW (that we do not live in a special place in the Universe), and apply it to a discretized 
matter distribution, then there is a way forward. In a seminal paper, Lindquist and Wheeler 
[55J showed that it is possible to construct an approximation to a cosmological space-time by 
considering a regular lattice of Schwarzschild space-times with metrics given by the line-element 



They did so for a closed universe and found that, in the limit of a large number of masses, such 
a space-time evolves almost identically to a closed FRW Universe: The separation between 
islands of mass increases, reaches a maximum, and then contracts with exactly the same time 
dependence as one would find for a perfectly smooth, dust filled universe. Configurations of 
gravitating discrete masses have also been considered in [56J. 

In this paper, we wish to explore the Lindquist- Wheeler model further, and to investigate 
the optical properties of such a universe. It is not the goal of this paper to propose a complete 
alternative to the currently favoured FRW universe dominated by a cosmological constant. 
Hence, we will restrict ourselves to considering only regular distributions of masses, and often 
only to the equivalent of an EdS universal. We expect that such an analysis should allow us 
to uncover some of the essential features of these space-times, even if it is not the most general 
case possible. We will return to the issue of extending our study to other spatial curvatures, 
matter contents and distributions of mass in future publications. 

The paper will proceed as follows: 

In Section [2] we recap the lattice model proposed by Lindquist and Wheeler in 1957 [55J. 
We reproduce their main result that the dynamics of a lattice in a closed 3-space obeys an 

3 Even though this assumes a matter content that is almost three times that which has apparently been ob- 
served through dynamical estimates with clusters and peculiar velocities, using FRW relations between redshift 
and distance. 
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evolution equation with the same functional form as the Friedmann equation, and with a scale 
that approaches the FRW one as the number of points in the lattice is increased We then extend 
Lindquist and Wheeler's study to include more general space-times, with arbitrary spatial 
curvature, and with a cosmological constant. We extend the 'cosmological time' coordinate 
they developed away from just the boundaries of the cells, giving us a global coordinate system 
with which to study the propagation of photons in, and between, cells. The problem of over- 
lapping Schwarzschild space at a boundary is discussed. The dynamics of the lattice is then 
shown to behave analogously to flat, closed, and open FRW universes, with or without a 
cosmological constant, when the curvature of the 3-space in which the lattice is constructed is 
chosen appropriately. 

With our approximate space-time in hand, we then proceed in Section [3] to investigate null 
trajectories within it. The relevant null geodesic equations are presented in the coordinate 
system found in Section [21 We explain how to propagate light rays between lattice cells using 
two different methods; one simple and approximate, the other more elaborate, and accurate. 
It is found that the results we present in later sections are reasonably insensitive to which 
particular method is chosen. 

We are then in a position to study redshift, as measured by observers looking along the null 
trajectories discussed above. In Section H] we present the formalism needed to determine these 
shifts, and make both analytic and numerical approximations to their solutions. The measured 
redshift to a radiating source is found, in general, to approach a value that is smaller than 
the corresponding shift in FRW, with relatively little scatter around the mean. The relation 
satisfied is 1 + z ~ (1 + zfrw) 7 ^ 10 - It is also found that deviations from this rule occur for 
trajectories that follow a special direction (such as being confined to a plane that picks out a 
symmetry of the lattice). 

In Section [5] we solve the Sachs optical equations along our null geodesies. These equations 
determine the expansion and shear of a bundle of null rays that are focused at either the source 
or observer. We provide a method of propagating these quantities between cells, and make 
analytic and numerical approximations of the distance measures that result. We also note the 
importance of the caustics that can develop due to the influence of shear. 

Section [6] contains a discussion of the observational consequences of the preceding sections. 
In terms of redshift, we find that the luminosity distance for a spatially flat dust dominated 
universe with A = takes the form tl oc (l + z) 2 — (l + z)~i , which corresponds to a deceleration 
parameter of go = 8/7 (i- e - objects at the same z appear brighter than in EdS). In terms of 
cosmological time, however, objects at the same t can appear dimmer. We present the Hubble 
diagram for this space-time, in the form of a plot of distance modulus, and continue to speculate 
on other cosmological observables such as CMB anisotropies, baryon acoustic oscillations and 
galaxy number counts. 

In Section [7] we conclude. A series of appendices then follow. 

Throughout the paper we attempt to present as many of our results as possible using analytic 
methods. We back this up with Monte Carlo simulations from a ray tracing code that allows 
us to confirm their accuracy. 
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2 An Approximate Space-Time 



Inspired by the success of the Wigner-Seitz construction in electro-magnetism [571EE], Lindquist 
and Wheeler (LW) constructed in [55] a lattice model of the Universe. The ideas they put 
forward are of central importance to our study, and so we reiterate them here. 

2.1 The Lindquist- Wheeler Model 

Starting with a positively curved hyper-sphere, LW distributed a number of "mass concentra- 
tions" into a regular lattice that they formed from tiling the 3-space with regular polyhedrsQ- 
Each cell of the lattice was given a central mass, and then approximated by a sphere, with the 
true geometry of the space-time (that which would result from solving Einstein's equations) 
being replaced by the Schwarzschild geometry of the closest mass. As LW noted: "This approx- 
imation demands that the distribution of gravitational influences just external to each sphere 
should depart relatively little from spherical symmetry" . 

The accuracy of this approach can be evaluated in the Wigner-Seitz construction by com- 
paring to known exact solutions [521 ED], with favourable results. Such solutions involve "empty 
lattices", where the potential is taken to be a constant throughout. In General Relativity a 
similar test can be performed, but this time with a cosmological constant dominating the grav- 
itational interaction |61j . The lattice results can again be seen to approach the exact solution 
(de Sitter space) in the appropriate limits, lending credance to the lattice model. 

The essential difference between the gravitational and electro-magnetic cases, apart from the 
non-linearity of the field equations in General Relativity, is that in the lattice model constructed 
by LW, the lattice itself is dynamical. This is because the non-zero normal derivative of the 
gravitational potential at the cell boundary results in a relative motion between the boundary 
and central mass. It is this motion that LW concerned themselves with, and with which they 
constructed an approximate global, dynamical space-time out of the Schwarzschild solution 
alone. 

For a cosmological interpretation, the scale factor of FRW cosmology now has to be replaced 
by some measure of the 'size' of the lattice, and the Hubble rate has to be replaced by the rate of 
increase in this size. To go further we must therefore explain what measure of size is intended, 
and with which time coordinate the expansion rate is defined with respect to. A key concept 
in the LW model, in this respect, is the idea of tangency between the constant curvature 
background hyper-sphere (on which the lattice is defined), and the 3-spaces that are locally 
orthogonal to the trajectory of the boundary sphere of each cell. This tangency gives a natural 
time coordinate, tlw, with which to describe the expansion (at least in the vicinity of the 
bounding spheres). What is more, it is a time coordinate that has some global meaning, as it 
can be used to define a congruence of time-like trajectories that are orthogonal to a common 3- 
space. As LW build their lattices in hyper-spherical 3-spaces, they also have a natural measure 
for their size: The radius of the hyper-sphere in an embedding Euclidean 4-space, cllw- 

It was found by LW that, in terms of the coordinates discussed above, the dynamics of their 
4 Tilings of 3-spaces of constant curvature are discussed in Appendix A. 
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Figure 1: The maximum radius of expansion, as a fraction of the FRW value, for lattices with 
N = 5,8, 16, 24, 120 and 600 cells. The dashed line corresponds to the maximum of expansion 
in a spatially closed, dust dominated FRW universe. 



lattice is specified by an evolution equation of the form 

2m 1 



a lw a lw sin3 V> a lw 



(4) 



where the over-dot represents a derivative with respect to tlw, m is the Schwarzschild mass 
at the centre of each cell, and ip is the (constant) angle subtended at the centre of the hyper- 
sphere between vectors in the Euclidean embedding 4-space that connect the centre of the 
hyper-sphere with the centre and spherical boundary of one cell. Clearly, this equation has 
the same functional form as the Friedmann equation, (j2J), with a dust-like energy content, and 
positive spatial curvature. 

The LW lattice therefore evolves in a similar way to a dust dominated, closed FRW universe. 
The only difference is the value of the maximum radius of the universe, when expansion ends 
and the onset of collapse is about to begin. In the LW case the maximum radius is given by 

(max) _ 2M 

a LW ~ ,r • 3 / ' v D / 

TV sin ip 

where M = Nm is the total mass of all cells in the lattice. In the FRW case we have 

(max) _ 4M , v 

a FRW ~ ' W 

where the total mass M is now written in terms of the energy density from the Friedmann 
equation as^| M = 2n 2 a^ RW p. The ratio of these two maximum radii is clearly independent of 
the total mass, M, in the space-time, and depends only on the number of cells in the lattice, 
N, and the angle they subtend at the centre of the hyper-sphere, ip. It provides a measure of 
the departure from FRW evolution, with the dynamics of the lattice approaching FRW in the 

i. (max) i (max) 1 

limit a> LW lo> Fm } -> 1. 



3 The surface area of a hyper-sphere of radius r is 2ir 



2^3 
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We know all possible values of N from Appendix A, and can straightforwardly work out the 
value of t/> for each lattice once it has been specified how the lattice cell polyhedra are to be 
replaced by spheres. The most successful case that LW consider is the direct generalization to 
curved space of the Wigner-Seitz approximation: That the boundary sphere of each cell should 
occupy the same fraction of the hyper-sphere as the original polyhedron^]. The results of this 
method are shown graphically in Figured! It can be seen that as N — > oo, and the continuum 
limit is approached, the lattice approaches the dynamical evolution of a spatially closed, and 
dust dominated universe. At N = 600 the difference is already less than 1.5%. This is the main 
result of LW. 



2.2 Geometric set-up 

The LW construction was engineered to allow the geometry inside individual cells to be approx- 
imated by the Schwarzschild solution, ([3]). Gluing these cells together in a suitable way then 
gives a global, dynamical lattice space-time. They found the dynamics of a lattice universe in 
a closed space, and showed them to be similar to the corresponding perfect fluid model. We 
now wish to generalise their model to include other spatial curvatures, and, crucially, to allow 
us to calculate the optical properties of the space-time. 

The Schwarzschild coordinates are, of course, perfectly acceptable when considering the 
geometry inside a single cell. They are not, however, well suited to describing the global 
geometry of a universe with many discrete sources. The problem is that these coordinates will 
not mesh at cell boundaries - that is, the 3-space of constant t from one cell will not be tangent 
to the corresponding 3-space from any other. This makes gluing cells together problematic, as 
the space that we would construct by putting cells at constant t next to each other would not 
be at all smooth: The coordinate patches for the spaces would intersect, rather than overlap, 
making any global interpretation of 'time' very difficult indeed. This problem is illustrated in 
Figure [2], with 3 intersecting 1-dimensional spaces placed next to each other. 

The problem was solved in LW by introducing a new time coordinate in the vicinity of the 
spherical cell boundaries that ensured the 3-spaces of constant time were orthogonal to the 
trajectory of the boundary. This allowed the two space-like surfaces of two adjoining cells to be 
tangent at any point where the boundaries met, hence allowing a common definition of time. 
This situation is illustrated in Figure [3J We now want to generalise their coordinates to spaces 
with different spatial curvatures, and away from the boundaries so that we can propagate null 
trajectories all the way through each cell in a consistent way. 

To do this, first consider a single Schwarzschild cell with a spherical boundary in free-fall. 
Now perform the following transformation from Schwarzschild time, t, to a new time coordinate, 
r: 

/ E - (1 - 2m) 

dr = v^Edt - V ^_2rn^ dr, (7) 

6 LW also consider a "less reasonable criterion" , in which the bounding sphere is taken to be just large enough 
to touch its nearest neighbours. This method was found to be considerably less successful than that mentioned 
above, and will not be considered further here. 



7 



Figure 2: A schematic of three one- dimensional space-like regions with intersecting, but not 
overlapping, surfaces of constant time, placed next to each other. At the boundary between cells 
different normal time-like vectors from different cells point in different directions. This is the 
way surfaces of constant t behave in our model. Hence, t does not represent a good choice for 
a global, cosmological time coordinate. 

time 



- space 

Figure 3: An example of three space-like regions with overlapping, rather than intersecting, 
surfaces of constant time. The normal time-like vectors from different cells can now be identified 
at the boundary. This is the way surfaces of constant r behave, and hence r can be interpreted 
globally. The smooth cell bottoms here represent the non-divergent behaviour of the r =constant 
surfaces at r = 2m. 




where E is a positive constant. The Schwarzschild line-element ([3]), describing the geometry 
inside the cell then becomes 



, n if 2m\ , „ 2 / / 2m\ , , dr 2 „ 
ds 2 = -- (l - —J dr 2 - -^E - (l - -jdrdr + — + rW. (8) 

In the limit E — > 1 this reduces to the well known Gullstrand-Painleve coordinates [62, 63J. 
The trajectory of a radially out-falling time-like geodesic is then given by 

*y =(*-i) +52, ( 9) 

dr J r 

where r is also the proper time measured along the trajector)Q. Each free-falling element of 
the boundary now has the 4-velocity 



u a = I l ; J(E-l) + — ,0,0 ) , (10) 



r 



and, for an arbitrary vector in the surface r =constant, given by n a = (0; n r , n e , n^) , it can be 
seen that 

u a n n = 0. (11) 



7 For in- falling trajectories one should take the opposite sign for the square roots in ([7]) and ([5]). 
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Figure 4: An example of the process of replacing regular polytopes with n-spheres of the same 
volume in a flat 2-space with square tiling. At the corners of the squares there is a "no man's 
land" that is not covered by the spheres, while at the centre of the edges of the squares there is 
an overlap between the spheres of neighbouring masses. 

The surfaces of constant r are therefore orthogonal to all infalling boundaries that satisfy (jUj). 
We will use r as our 'cosmological tim^fl'. 

We can now consider how to construct a lattice from individual cells with the geometry 
specified by (jHJ), and spherical boundaries at r = a(r) that satisfy (jHJ). To do this we must first 
choose one of the regular lattices specified in Tabled] of Appendix A, at some initial time tj. A 
mass, m, must then be placed at the centre of each cell. This prescription leads to the highly 
symmetric lattice structures that we wish to study her^|. 

We now want to make a Wigner-Seitz-like approximation, following LW, and replace the 
polyhedron lattice cells with spheres. Generalising 'Condition I' of [55] , we choose to specify 
these spheres as occupying the same amount of spatial volume as the original polyhedra in the 
constant curvature space in which the lattice is defined^!. This will require the spheres from 
neighbouring cells to overlap at the centre of the faces of the original lattice cells, while leaving 
a "no man's land" at their corners (a region external to all spheres). This type of situation is 
illustrated using the square tiling of a flat 2-space in Figure SJ 

8 We use the letter r, as this is proper time of observers following trajectories given by (|10p. It is not to be 
confused with the 'conformal time' coordinate often used in cosmology. 

9 One could be less prescriptive about the symmetries involved in tessellating the 3-space, leading to more 
general situations. We will discuss this briefly later on, but will postpone a detailed investigation for future 
studies. 

10 This is equivalent to the Condition I of [55] on lattices constructed in positively curved 3-spaces: They will 
occupy 1/N of the total solid angle on the hypersphere. 
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If we now replace the geometry inside each cell by (jHJ) then the spheres will have a motion 
specified by (GO), and the orthogonality requirements we considered above will lead to tangency 
of space-like surfaces of constant r at the points where the spheres intersect. Away from these 
points, in the 'overlap' and 'no man's land' regions, this tangency will no longer exist. It will, 
however, be broken in an opposite sense in each of the two types of regions, and so will, in some 
sense, be satisfied 'on average'. 

The dynamical properties of these space-times can now be determined. 

The case of a space-time containing a cosmological constant is investigated in Appendix B. 

2.3 Cosmological evolution 

From Equation it can be seen that for E < 1 the boundary spheres will reach a maximum 
of expansion, at a max = 2m/(l — E), and then begin to recollapse. For E > 1, however, the 
spheres will expand eternally, and reach infinity with a velocity da/dr > 0. The dividing case, 
E = 1, corresponds to shells that just reach infinity, and hence to the escape velocity. Clearly 
this behaviour is analogous to the behaviour of FRW cosmology, with closed, open, and flat 
spatial curvatures, respectively. This analogy is made explicit if one re-labels variables such 
that m —>■ M/2 and E — > 1 — k. In this case ([9]) becomes 

4 = 4-4. m 

a 1 a 6 a 1 

which is obviously just the Friedmann equation, ([2]), with p = 3M/8na 3 . Henceforth, we will 
use a to label the position of the boundary sphere, and r as our radial coordinate inside each 
cell. The volume of the lattice cells at future times is then given by the corresponding volume 
of the bounding spheres. 

We can make contact with the LW treatment by noting that in the case of a closed space-time 
the LW variable ip can be related to E via E = cos 2 ip. Clearly, cos 2 ?/' < 1, and so they found 
their space to reach a maximum of expansion, at the value found above, before recollapsing. 
The generalisation of this to hyperbolic space is then E = cosh 2 ip, which is always > 1, leading 
to the hyperbolic expansion found above. 

As in the study of LW, we find that all of the models we consider are governed by evolution 
equations with identical functional form to FRW space-times. It was then shown in [55J that for 
the case of spatially closed universes, while the evolution equation had the same form, the scale 
of the problem did not: The maximum of expansion was different in the two cases. Here we will 
focus on the case of a spatially flat universe, with E = 1. The expansion is then scale-invariant, 
and so we do not anticipate any such discrepancies between the discrete and continuous cases 
to occur. 

Now, while the LW construction only considered surfaces orthogonal to the trajectory of the 
boundary, in the neighbourhood of the bounding sphere, the space-like surfaces r =constant, 
in the coordinate system (jSj), are orthogonal to all shells in radial free-fall that obey This 
allows us to extend the meshing coordinates away from just the boundary region, and gives 
us a time coordinate that we can consider to be the analogue of cosmological time in perfect 
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fluid cosmolog y. The idea of the single bounding sphere in LW is here replaced by a continuous 
family of shellg^j following a congruence of time-like geodesies specified by ( fTUl) . and which are 
identical to each other, up to a translation in r. Such a coordinate system will be very useful in 
determining how time-like geodesies should pass between cells. This is particularly true in the 
spatially flat universe, in which case the surfaces of constant r, in (jHJ), are simply the Euclidean 
3-spaces in which the lattice cells are defined^!. We will concentrate on this case in the analysis 
we perform in later sections. 



2.4 Validity of the model 

The model we are considering, as outlined above, has features that appear at first glance to 
be simultaneously more realistic and less realistic than standard FRW cosmology. The raison 
d'etre for this model is that it does not a priori assume that discrete objects can be simply 
approximated by a continuous energy density; it explicitly maintains the discrete nature of the 
objects, and contains the continuous approximation as a limit. In this regard we consider it to 
be a considerable improvement on the more usual perfect fluid description. 

Having said this, there are also obvious draw-backs. The model itself relies on certain 
approximations, such as the space-time being only an approximate solution of Einstein's equa- 
tions. In the electro- magnetic analogy these approximations have been shown to be well justified 
[591 160] . Similar tests suggest this is also true in the general relativistic case [6TJ , though further 
study is needed to determine this conclusively. 

Beyond this, one may also question the validity of approximating the matter content of 
the Universe as identical spherically symmetric mass distributions that are equally spaced on 
a regular lattice. At first glance this is clearly not true of our Universe. Whatever unit of 
structure we conceive of, it certainly is not arranged on a regular lattice. Detailed surveys of 
the distribution of galaxies have revealed an intricate tapestry of nodes, clusters, filaments and 
walls that seem to have resulted from a stochastic process of structure formation. In fact, it 
has been shown that the morphology of these structures is fractal over a wide range of scales 
[65j-[68j. This suggests that a less symmetric structure than a regular lattice would be more 
realistic. However, we wish to retain as much as possible of the cosmological principle^! in our 
study, and it makes sense to consider the simplest, most symmetric case first. 

Now let us return to the question of the masses. Again, most objects in the Universe do 
not appear to be identical and perfectly spherical. Never the less, it is certainly not unusual 
to approximate the gravitational fields of non-spherical objects (such as disc galaxies) as being 
spherical, and taking equal mass objects seems like a justifiable approximation in order to make 
progress in understanding the problem at hand. Taking into account the detailed shape and 
structure of every object in the Universe would obviously be prohibitively difficult. One would 

11 We still keep the concept of a single 'special' shell at r = a that we will use to determine the volume of the 
cell. 

12 For the closed universe models, with E < 1, the coordinate system given in ([8|) has the problem that it 
does not cover all of the 'no man's lands' in the vicinity of the maximum of expansion. In this case a more 
complicated coordinate system can be found, as shown in [55] , 

13 A statistical cosmological principle (in which one assumes statistical homogeneity and isotropy) would be 
more desirable, but is currently impossible to implement. The next best thing appears to be a regular lattice. 
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hope that future studies would make progress from the simplest case studied here, to more 
realistic situations. 

An obvious choice for the central masses in this model is to consider them to be galaxies 
with the same mass as the Milky Way. In this case we find the spacing between masses should 
be of the order of a few mega-parsecs. Clusters of galaxies would be another sensible choice, 
and would correspond to a fractionally larger spacing of masses. Alternatively, one might also 
consider dark matter particles (with masses of a few 100 GeV). Such a choice may be more 
suitable for applying the model considered here to the early Universe. The lattice spacing would 
then reduce to a few metres. 

Throughout this paper we will use Milky Way type masses, although we have also considered 
other objects and found our results to be largely insensitive to this choice (unless we consider 
very large masses). 



3 Photon Trajectories in a Lattice 

We will now consider the trajectories of null particles in this space-time. These will be of basic 
importance for understanding the optical properties of the lattice universe. All analytic results 
presented here will be backed up by numerical analysis. 



3.1 Null geodesic equations 

The Euler-Lagrange equations derived from the line-element (jHJ) are 



£ (l-^Ws-Cl-^Vl = (13) 



dX \ \ r J W \ r J 



dX 

together with the null constraint 



dr / _ / 2m \ df m 9 „a 9 ,_, 9 „ ; 9 , ,s 

— -JE-II — = — -f 2 + Er9 2 + Er sin 2 9(f) 2 14 

dX W \ r J d\ r A 

4r(r 2 e] = t 2 sin 9 cos 9(f) 2 (15) 
dX V / 

11 V sin 2 00) = 0, (16) 



_ ^i_^f 2 + r 2 + Er 2 9 2 + Er 2 sin 2 9<j> 2 - E- (l-^jrf = 0, (17) 

where dots indicate derivatives with respect to A, an affine parameter along the geodesic. 
Integrating these equations gives 
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B 

J<b 



(19) 



r 2 sin 2 9 

' (20) 



r 4 r 4 sin 2 
5 2 J 2 / 2m 



where -B, J and are constants. Clearly, one could rotate coordinates so that 9 = n/2 in each 
cell, although this will not be as useful here as it usually is, as we will want to match different 
coordinate systems between cells. Rotating them differently in each cell would confuse things. 

The geodesic equations can also be found in Cartesian coordinates (when E — 1), which 
are useful for numerical implementation, and in the presence of a cosmo logical constant. The 
equations for these cases are given in Appendices C and D, respectively. 



3.2 Matching trajectories at boundaries 

Now consider the following situation: A photon is emitted in one lattice cell and passes through 
a cell boundary (or, more likely, several cell boundaries), before it is observed. In order to make 
predictions about observations of events that occur in cells that are at some distance from the 
cell inhabited by the observer, we must be able to propagate photon trajectories between cells. 

One may initially suspect that at the boundary between cells it may be suitable to simply 
perform a transformation of spatial coordinates from one cell to another, via a suitable transla- 
tion of the origin of the coordinate system. However, it is soon seen that this is not a viable way 
of propagating null trajectories between cells. To see why consider the null constraint equation, 
(I17p . and a photon that is following a radial trajectory in both the first and second cells. If we 
denote the spatial coordinates of the first cell by un-hatted coordinates, then ffT7|) gives 

- (l - H=) * + * - 2^-(l-f^ = 0, 

where r > 0, as the photon is leaving the cell. If one were to perform the simple translation 
discussed above, keeping f the same, then the equation above would transform to 

-(l-£)f + f" + 2^-(l-£)* = 0, 

as r = — f and r = f at the boundary. Hats denote the spatial coordinates of the second cell. 
It can immediately be seen that this new equation is not compatible with the null constraint, 
(|17p . in the second cell, as it has the wrong sign before the third term on the left-hand side. 
We must therefore be more careful. 
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Our basic criterion for matching trajectories across boundaries is that observable quantities, 
such as photon frequencies and directions, should be independent of which coordinate system 
an observer moving with the boundary chooses to use. 

Now, at the boundary, it will only be a set of measure zero trajectories that actually pass 
through the single shell that is described by our 'special' bounding sphere, which prescribes the 
volume of the cell. All other trajectories will either pass though in the 'overlap' region, which is 
the intersection of the interiors of the two bounding spheres, or in the 'no man's land', that is 
outside of all bounding spheres. This means that observers comoving with the family of shells 
that obey ([TO]) will, in general, be in relative motion with respect to the cell boundary when 
the photon in question passes them. Such motion should be expected to result in a redshift 
between comoving observers from neighbouring cells who are both at the same point on the cell 
boundary at the same time. We will now consider two different methods of accounting for this 
effect. 

In Method I we will appeal to the 'average' tangency between the space-like volumes of 
constant r. A simple, approximate matching criterion is then given by the condition that f is 
the same on leaving one cell, as it is on entering the next. This will be shown in subsequent 
section to correspond to the condition that, at the boundary, the frequency of a photon measured 
by a comoving observer from the first cell should be identified with the frequency measured by 
comoving observers from the second. Clearly, this criterion will not be exactly satisfied by all 



through the boundary in the 'overlap' region will pick up a change in f with the opposite sign 
to those that pass through in the 'no man's land'. 

In Method II we will attempt to account for the relative motion between comoving observers 
from different cells in a more detailed way. Instead of assuming that the effects described 
above cancel on average, we will work out the redshift between a comoving observer who is 
momentarily at the cell boundary, and an observer who is moving non-geodesically with the 
boundary. Rather than identifying the frequency measured by comoving observers from each 
cell, we will then identify the redshift measured by the observers who are moving with the 
boundary. This method will be more complicated, but we expect it to account for the imperfect 
tangency in a more complete way. 

The results of using Method II will turn out to be very similar to the results of using Method 
I. This provides us with motivation for considering the concept of average tangency between 
the space-like surfaces of neighbouring cells to be a valid one. 

3.2.1 Method I: A simple, approximate matching 

Using the matching condition that f on leaving the first cell is the same as f on entering the 
second cell, we are left with the task of finding f, 9 and <f> in the new cell. To do this consider 
the following decomposition of the 4-vector tangent to the null geodesic, 




|, but it will be satisfied 



in an average way over many trajectories: Those that pass 




(22) 



In fact, any that do not pass through the 'special' bounding sphere at the boundary. 
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where u a is the trajectory of an observer on one of the free falling shells specified by (II 01) . The 
4- vector n a is the same unit space-like vector that was considered below Equation (1101) . 

Now consider describing the null geodesic in a neighbourhood of the observer using the 
coordinate system of the first cell. This results in 



k a = t\l-\l{E-l) + ^ + n\n\n*\, (23) 

where we have made use of the fact that u a k a = — f , and (n r ,n e ,n^) obey the normalisation 
condition 



+ r (n ) + r sin 9{rv) — 1- (24) 

Now consider the same geodesic, in the same neighbourhood of the same observer, but this 
time using the spatial coordinate system of the second cell. In this case the same reasoning 
gives the same tangent vector as 



k a = f\ l;^(£-l) + y+n f ,nV* | , (25) 

where hats denote indices in the second coordinate system (that of the new cell into which the 
photon is propagating). As discussed above, the cell boundary is equidistant to the central 
mass of each cell, so that r = f, and the vector n a obeys a similar normalisation condition to 
that of n a . 

Method I now tells us that f in the neighbourhood of the observer with 4-velocity (11 01) at 
the boundary is, on average, common to both coordinate systems, so that 

r|out - r|in. (26) 

We can then relate , 9, (pj to (r, 9, (f\ using the condition that the projection of the tan- 
gent vector of the null geodesic into the rest space of an observer at the boundary should be 
independent of which spatial coordinate system is being used. 

To see how this works consider first a radial geodesic. The condition 9 = <p = then 
gives n e = r& = 0, and similarly in the hatted coordinates. We also have that n r = 1 in the 
coordinates of the first cell (this says that the geodesic is outgoing). In the coordinates of the 
second cell, however, we have n r = —1, as the geodesic is ingoing in the new coordinate system 
(i.e. moving to smaller r). We can then write the relation between r\- m and f| out as 



2f ] /(E-l) + 2 ^. (27) 

This satisfies the constraint equation exactly, and we avoid the problem noted above. 

Now consider a general (not necessarily radial) geodesic. In this case we have that the 
expression for k a above gives, for E — 1, 

9 = fn 9 (28) 

= fn^ (29) 

/2m , , 

r = fn r + fd—. (30) 
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Similar relations are obeyed by the hatted coordinates. The expression k a k a = now gives 

k a k a = k a k b g ab = (n a n b g ab - l)f 2 = {n a n"5 a p - l)f 2 , (31) 

where Greek indices run over spatial coordinates, and 8 a p represents the spatial metric of 
the Euclidean 3-space. This clearly satisfies the null condition, as n a n a = 1. The constraint 
equation can then be seen to be valid in any two coordinates systems related by a transformation 
of spatial coordinates of the form x a = A a a x a , where A 7 Q,A a a = <5 7 a . 

The procedure for propagating a general geodesic through the boundary using Method I is 
then the following: 

(1) Take t,f,9 and <fi from the first cell and use them to find n r ,n e and n^, using the 
relations (12gj)-(13TJj). 

(2) Transform the (r, 9, 0) coordinate system into a system with the same symmetries as 
the cell (this is Cartesian coordinates for E = 1). Translate the origin of the coordinates from 
the central mass of one cell, to the central mass of another. Transform to a new set of spherical 

coordinates, ^f, 9, 0^ . Calculate n r , n® and in these new coordinates. The vector n a is then 
given in terms of n a by the transformation 

" a = a?"* (32) 

where the two coordinate systems are related, for E — 1, by 

f 2 = r 2 + Xq — 2rxo cos cp sin 9 (33) 
2 ? (rcos0sin# - x ) 2 , . 

C ° S = (x 2 + sin 2 9-2r Xo sin 9 cos 0) (34) 

2 9 C ° s2 ^ (35) 

(r 2 + Xq — 2rxo sin 9 cos 0) ' 

and Xq is given by the translation between Cartesian coordinate systems, x = x — Xq. This will 
be equal to the width of a cell, at the time the photon hits the boundary. 

(3) Use these values to calculate f, 9 and <fi using fl28p - (|3H|) . given that f = f. These are 
then the initial conditions for propagating the null geodesic through the new cell. 

(4) Repeat at the next cell boundary. 



3.2.2 Method II: A more elaborate, and accurate matching 

The above method did not take into account the fact that along individual trajectories the 
observers with 4-velocity u a on either side of a cell boundary can have a relative velocity 
between them. This is due to the imperfect tangency between surfaces of constant r on either 
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side of a cell boundary. To calculate the effect of this we will consider an observer moving 
(non-geodesically) with the boundary. The 4-velocity of such an observer, w a , can be given by 

w a = (w T ;w r , 0,0), (36) 

and satisfies w a w a = — 1, if the observer is time-like, and radially moving away from the mass 
at the centre of each cell. The normalisation condition gives us that 



^ - ' • (37) 

The 3-velocity, w r , can then be worked out straight-forwardly. For E = 1 we can use Euclidean 
geometry to give w r in terms of the expansion of the shell that defines the volume of the lattice 
cell, d, as 



"K 

" ; 6 



1/3 



\Jl + (5x) 2 a, (3* 



where 5x is the distance from the centre of the face of the particular lattice cell that is being 
crossed, as a fraction of the distance from the cell centre to the centre of the cell face. 



The frequency of a photon with tangent 4-velocity k a measured by this observer is then 
given by 



- w a k a = f\l 1 - — + (w^ v 2m , r. (39) 



- + {wr)2 -ihm 

The values of f and r are given by the solutions to the null geodesic equations stated above. 
Similarly, in the second cell we have 



■W- ( w r _ /am X _ 2m + / r \2^ 

2 ™ + f w r)2 _ 1 V ' V L 

r 1 j (1-^) 



w a h =r X /l-— + (w r ) 2 v ji 2m , f, (40) 



where, by symmetry, the 4-velocity of the observer in the cell boundary is the same for both 
cells. We now wish to identify the frequency of the photons measured by the observers moving 
with the boundary, giving the condition —w a k a = —w a ka- 

If, as before, we now use the projection of the 4- vector k a into the Euclidean 3-space, then we 
have that 9, <j) and r in the first cell are again given by Equations fT28l) - fl30l) . Similar expressions 

are satisfied by 9, <j> and f, but with hatted coordinates. 

We can now use fl3"9l and (|40l) to calculate f in terms of known quantities, as 



* = V 7zr, ; "\ t (41) 



For w r = \j2m/r we recover f = f, as used in Method I. More generally, the value of f given 

by (14"T|) can be used, together with n 6 , r& and n r , to find 9, and f. The steps to be followed 
at the boundary using Method II are therefore: 
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(1) Find vf at the boundary, using ( 1381) . 

(2) Take t,f,6 and <fi from the first cell and use them to find n r ,n and n^, using the 
relations (|28|)-(|30|). 

(3) Find n a in terms of n a via the transformation 

n> = 5F»*- < 42 > 

where the two coordinate systems are related by the same coordinate transformation as in 
Method I. 

(4) Find f using (Hil . 

(5) Use the values found in (3) and (4) to calculate f , 6 and using (|28l) -( |30l) (with hats 
added). These are then the initial conditions for propagating the null geodesic through the new 
cell. 

(6) Repeat at the next cell boundary. 



3.3 Numerical Implementation 

We have implemented the methods described above numerically, so that we can propagate 
geodesies with any initial positions and directions out to arbitrarily large distances. This is 
effectively the equivalent of the ray tracing methods applied in simulations of lensing through 
cosmological configurations. Note that in all plots, we integrate back in time from the observer's 
current position, and fix the affine parameter, A, to be initially zero. 

A few brief comments on this seem in order at this point. We find that while integrating 
the geodesic equations through a Schwarzschild cell it is preferable to use Cartesian coordinates 
(as given in Appendix C), rather than the natural, spherical coordinates of the Schwarzschild 
geometry. In this way we can avoid coordinate singularities (where, for example, sin# = 0). 
Although these singularities are rare for any single geodesic, they can still occur occasionally if 
one integrates over many different geodesies in order to obtain statistically significant results 
from a Monte Carlo simulation. 

A further point to comment on is that with our choice of Milky Way sized masses, the 
ratio between the Schwarzschild radius and the cell size is minute (of order 10~ 8 ). This means 
that large deflection events (as trajectories pass nearby the central mass) are almost non- 
existent, and that all but a few geodesies suffer only negligible deflections as they traverse a 
single Schwarzschild cell. With this fact in hand, we find it sufficient to consider rectilinear 
trajectories (akin to the Born approximation), and to calculate our optical quantities along 
them. 
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4 Cosmological Redshift From Schwarzschild Patches 



Now that we have the equations for null geodesies, we can calculate the redshifts between source 
and observer that are so important in cosmology. The redshift will, of course, depend on the 
motion of the source and observer. These can, in principle, be completely arbitrary. However, 
the closest analogy to a comoving source and observer in FRW cosmology will be a source and 
observer that are comoving with one of the family of free-falling shells with 4-velocity specified 
by (fIQ]h As always, the redshift, 1 + z, is given as the ratio of the frequency, —u a k a , measured 
at the source and observei@. Here this is 

- u a k a = r (43) 

so the redshift is 

l + z=?K (44) 
r\o 

where subscript e and o denote when the photon was emitted and observed, respectively. In 
general this quantity will need to be calculated numerically. We find, however, that we are able 
to deduce reasonably good analytic approximations to the numerical results. We will present 
our analytic approximations first, and then proceed to compare them to numerical solutions. 



4.1 Analytic approximation 

To find analytic approximations for the redshift consider first a single null trajectory. For an 
expanding space-time we then have from the geodesic equations that 



(l - a^f) 



r = B V n lm, ' ■ (45) 



The factor of a is included to account for the magnitude of r as a fraction of the total 'velocity' 



r 2 + r 2 9 2 + r 2 sin 2 6(j) 2 ~ B, and is given by a = r/B for any particular cell 16 !. 
For 2m <rwe then have 

f ~ B(l - v {a-v)), (46) 



where v = ^J2m/r is the velocity of the observer. The redshift induced by a photon travelling 
through a single cell, (1 + Szi), is then 

(1 + SZi) = ^ * 0- ~ ViniOin - V in )) „ 1 + _ 
T | out (1 — V out (a ou t - V out ) 

as B is constant inside each cell. In the last equality we have used V{ n ~ v ou t for a photon 
passing through a single cell. For small z the total redshift induced by travelling through n 



15 This is the case for both Methods I and II. The 4-velocity of the non-geodesic observers in Method II, w a , 
is only used for propagating photons across boundaries. 

16 For a geodesic with J = = we have a = 1, as the photon is moving entirely in the radial direction. 
More generally a G [0,1]. 
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Figure 5: The approximate relation expected between redshift in a lattice universe, z, and the 
corresponding redshift in an FRW universe, zp^w ■ The solid line is for null trajectories not 
lying on a principle axis, and the dotted line is for near radial geodesies. 



cells is now 



l + z = f[(l + Sz t )^l + 2j2^J^^l + 2J^T 1 /3 J (48) 

i=l i "' 

where denotes a at the moment the photon enters the ith cell. 

We now need i = i{r) in order to proceed (this is the cell number that the photon is in, as 
a function ofr). Assumed the expansion of the cell is slow compared to the scale of the photon 
crossing time 17 !, to lowest order in m/r we then have that the crossing time of a region, Ai, is 
given by 

^ -2a,, (49) 

with the continuum limit 

d?-- 2o »y ■ (5o) 

More generally, for non-radial geodesies, we can write 

dr f r \ 2/3 

— ~ -2&ao - , (51) 



di \tq J 

where j3i is the distance across the cell as a fraction of the length of the radial trajectory. 
Substituting this into (14"8"1) and integrating gives the redshift as 
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This should be a good approximation as long as the cells are small. 
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Figure 6: The relation between redshift in a lattice universe, z, and the corresponding redshift 
in an FRW universe, zpjmr, for 5 different trajectories starting off with different orientations. 
Close to the observer {left panel) there is still a fair amount of scatter, but after traversing a few 
Schwarzschild domains, this scatter becomes negligible {right panel). We also plot the analytic 
fit with (7) = 2/3 from expression (f56j) {dashed line). Equation ff5Z>l) with (7) = 7/10 {solid 
line) lies directly on top of the numerical results {dotted line). 



V2m / a \ . r e . . 

1+2 - '-^{nr* (52) 

(t \ 2/3 

l-<7>lnpj (53) 

~ i + ( 7 )i n (^M (54) 
\a e J 

= 1 + (7) In (1 + z FRW ) (55) 

^ {l + z FRW )^\ (56) 

where (7) = (a/j3). Strictly speaking, the last step here is only valid up to linear order 
in In (1 + zfrw)- This is perfectly acceptable, however, as the calculation before hand has 
only been performed up to this order of accuracy. Now, if we use the fact that redshifts 
should combine as factor^, then we have good reason to suspect that ( 1561) will be a better 
approximation to 1 + z than (|55l) . The reason for this is that Equation (|56|) can be seen to 
combine factors of redshift in exactly the way required, and to have the correct linear term 
in an expansion of In (1 + Zfrw) around 0. For these reasons we expect it to be a reasonable 
non- linear extension of the relation between z and zfrw- This is confirmed numerically 
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Such that (1 + zi^ 3 ) = (1 + Zi-> 2 )(1 + ^2^3)- 
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Figure 7: The error in the z — zfrw relation for 40 random trajectories. 

The angle brackets here refer to an ensemble average, over all the cells along a trajectory. 
For radial geodesies (7) — > 1, and so z — > zfrw- F° r general geodesies it is shown in Appendix 
E that (7) ~ 2/3. The redshift in a lattice universe is then found to go like 

1 + z ~ (1 + z FRW ) 2 '\ (57) 

for a trajectory that is not aligned with any principle axis. This result marks a significant 
deviation from the corresponding observable in FRW. It is shown graphically in Figure [51 

4.2 Numerical Results 

If, instead of making approximations, we integrate the geodesic equations numerically, then we 
can get more accurate results, and determine the validity of the simple expression ( |57j) . Doing 
this, we find that the analytical results are good (though not perfect) approximations. As one 
would expect, we do find an initial scatter in the z — zfrw relation as the geodesies traverse 
the first few Schwarzschild cells. 

In Figure [6] we plot the relationship between z and zfrw for five trajectories with different 
initial directions. The steepest curve in the left of these plots is along a principal axis, and is 
clearly different from the others, with z ~ zprw- The directions of the other trajectories have a 
range of random angles and, although they start off with different gradients, they very rapidly 
converge upon the mean. As shown above, deviations from the mean depend on the number 
of domains traversed. In Figure [7] we plot the relative deviations from the mean, and find that 
for Milky Way sized masses they do indeed become negligible very rapidly. There is very little 
scatter at z > 0.1. 

We can now use our numerical analysis to show the limitations of the analytic approxi- 
mations. These are shown for comparison in the right-hand panel of Figure O Although not 
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Figure 8: The relative difference in the z — zfrw relation for the two choices of matching 
conditions at the cell boundaries. 

perfect, it is can be seen that (1561) with (7) = 2/3 offers a reasonable approximation. With 
(7) = 7/10, however, we find that Equation (1561) fits the numerical results almost perfectly. We 
will therefore use this semi-analytic value for (7) in our calculations, later on. 

It is also now possible to compare the two matching schemes we outlined above, for propa- 
gating geodesies between cells. As advertised, we find that the effect on physical quantities is 
negligible, showing that the 'average' tangency of 3-spaces is indeed a valid concept. In Figured] 
we show this explicitly with a plot of Az/z as a function of zfrw, f° r the two different choices. 
The relative error is at the level of less than 1% out as far as Zfrw — 2, which is remarkably 
accurate. 

In what follows it will also be useful to know the cosmological time, r, and redshift, z, 
as functions of the affine distance along the geodesies, A. One finds, both analytically and 
numerically, that 



where r is the current age of the Universe and A has been chosen so that it is at r = r , 
and — 1 at r = 0. In Figure [9] we plot both r = r(A) and z = z(X) all the way from the initial 
singularity (A = — 1) to the present day (A = 0). These results are derived from averaging over 
a set of 40 realisations. We find that the numerical results fit Equation f[5"8"l) almost perfectly. 



T 



(1_|_ \)2H+S , 



(58) 
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Figure 9: The functional dependence for both z = z(X) and r = t(X), averaged over a set of 40 
runs. 

5 Optics of a Discretized Matter Distribution 

While the dynamical evolution of a cosmological model is important for understanding, for 
example, the growth of structure, one could argue that it is of equal, or even greater importance, 
to understand its optical properties. While the latter requires a knowledge of the former, there 
is no guarantee that similar 'average' dynamical properties in two space-times should lead 
to similar optical properties. In fact, we have already seen that the lattice space-time gives 
different redshifts to FRW, and, as we discussed in the introduction, there are good reason to 
suspect that there should be considerable differences in distance measures too. 

There are various reasons to suspect this. Firstly, the dynamical equation for the expansion 
of a bundle of null geodesies is driven by a term oc R, the Ricci curvature tensor. In the 
continuous fluid approximation this is non-zero, while in a universe with discrete islands of 
matter it will be exactly zero everywhere outside of the matter itself. What is more, the shear 
in a homogeneous and isotropic perfect fluid universe is zero, while in a universe with discrete 
matter sources it will be non-zero. These differences have competing effects on the luminosity 
of astrophysical sources, the former making them dimmer in the lattice universe, and the latter 
making them brighter. These phenomena have been described by Bertotti [69J. 




To investigate the optical properties of the space-time we will have to integrate the Sachs 
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optical equations [51] along each null trajectory. These equations read 



^- + 9 2 -LU 2 + a*a = -]-R ab k a k h (59) 
dX 2 

^- + 2u9 = (60) 
^ + 2a6 = C abcd (t*) a k b (t*) c k d , (61) 

where 9, uj and a are the expansion, rotation and complex shear scalars, respectively. The C a bcd 
is Weyl's tensor, R ab is the Ricci tensor, and t a is a vector that is orthogonal to k a , null, and 
has a magnitude of 1 (i.e. t a k a = 0, t a t a = and t a (t*) a = 1). In Schwarzschild space-time 
R a b = 0, and the driving term in the a equation can then be found to be 

3m T 2 

C = C abcd (t*) a k h (t*) c k d = -jf-e w , (62) 

where ^ is a constant, specifying the complex phase. Once the expansion scalar is known, then 
the angular diameter distance is given by the integral 



ta oc exp | J 6d\ j , (63) 

and the luminosity distance is given by Etherington's theorem [70] as 

r L = (1 + z) 2 r A . (64) 

It can also be seen that for a point source u = is always an integral of (JBTij) . Choosing units 
appropriately we then have the two equations 

1 d 2 rA * , N 

— - rrf + aV = 65 
ta d\ l 

do 2adr A 3mJ 2 ^ 

rfA r A dX r 5 v 7 

The initial conditions for integrating these equations are then a\ Q = 0, r A \o = and <ir^/(iA| = 
constant. (The reader should not confuse r and r A in these equations). The value of r = r(A) 
should be substituted from the solution to the geodesic equations, together with the relevant 
value of J and \I/ in each cell. 

Let us now write the complex shear as 

a = |a|e i$ , (67) 

and define the new variable X = |c|r^. This expression for a has two independent real parts, 
\a\ and $, and therefore contains all the information about the shear scalar. 



In terms of these new variables the evolution equations for the optical scalars, (|65|) and (166 
become 



d 2 r A 
d\ 2 



+ X 2 = 0, (6* 
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dX 


3mJ 2 r\ 


~dX ~ 






3mJ 2 r\ 


X d\~ 





and 

AY ZrnJ 2 r 2 . 

cos(^ - $) (69) 
sin(^-$), (70) 

for the real and imaginary parts of ( |66i) . respectively. These last two equations can then be 
integrated to give 

X sin(^ - go) 

X sin(tf-$) ' 1 j 

where X , \l/o an d $o are constants. This results in the single equation for the evolution of the 
shear 

(72) 



fdXV = ( 3mJ 2 r 2 A Y f _ X 2 sm 2 ^ -$ ) 
V dX J v r5 / V ^ 2 



where the sign of dX/dX is given by the sign of cos(\& — $). Once we know Xq, $o and \l/o in 
a given cell, we can then solve the equations (j6"81) and (1721 . above. 



5.1 Evolving shear between cells 

We now want to know how to evolve shear between cells. It should be the case that observable 
quantities are continuous along trajectories, so we know that and X should be the same on 
starting the next cell as they were on leaving the last one. It remains to determine how the 
phase factors $ and \1/ are propagated across cell boundaries. 

It is clear from (|68p - (|72|) that the difference between $ and \I/ is of primary importance for 
evolving the optical scalars. To make clear what this quantity corresponds to, let us consider 
further the 4- vectors t a . We can think of them as giving a set of two mutually orthogonal 
space-like unit vectors, that are orthogonal to k a , via 

t a = ^=(p a + tq a ), (73) 

such that p a p a = q a c[a — 1 and p a k a = q a k a = p a q a = 0. The phases $ and ^ therefore 
correspond to the orientation of the shear and of the driving term (162]) in the plane spanned 
by p a and q a . The quantity ^ — $ is their relative orientation. 

Now, in a space-time that is globally Schwarzschild (i.e. with only one point-like mass in 
the whole space-time) it can be seen from (1701) that $ is driven toward This means that 
the direction that the beam is sheared in is drawn into alignment with the term driving the 
shearing, just as should be expected. In the present situation, however, when we pass between 
cells the orientation of the shear that has accumulated up until that point, $, will not, in 
general, be in the same direction as the term driving the shear in the new cell, \& '. Stated 
another way, the beam is sheared in different directions as it passes through the different cells. 
This change of direction is given by the change in \& — $ as the beam passes a cell boundary. 

We are, of course, free to perform an arbitrary change of t a as we move between cells. If we 
were to do this, however, then we should expect to have to work out both $ and ^ with respect 
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to the new set of 4-vectors. Instead, it makes sense to not rotate p a and q a at the boundary, 
so that we can keep the phase of accumulated shear up until that point, $, the same. This 
can be achieved by a direct transform from the coordinate system of the first cell, x a , into the 
coordinate system of the second, x a , via x a = x a (x a ). The vector t a then transforms as 

(1 1™ ^ 

f&\ _ _ + h \ — + a | (7A\ 

| entry to new cell Q^> l ex ^ f rom °ld cell — L |exit from old cell- V ' 

This is the same coordinate transformation that was used previously on n a , and hats mean 
quantities in the cell with the hatted coordinate system. We now need to know the phase of 
the driving term in the second cell, with respect to these vectors. 

To find this, consider that for the Schwarzschild geometry ([8]) and tangent vector (1251) we 
can write t a , for E = 1, as 

t a = e-^ /2 {t a + ak a ) } (75) 



where 



and where \1/ and a are real and complex functions, respectively 19 !. To see that this \1/ is the 
same as that in (I62p . simply contract it with the Weyl tensor using the result@ 




t" ~ T \ irf; + \l—n r \ , - sin Qn\ ) . ( 76 ) 



C = C abcd k a (i*) b k c (i*) d = -Lf-. (77) 



We now need to know how \I/ and a evolve inside each cell. To find this let us impose that t a 
is parallelly propagated along the null curve. This then gives 



d\ ; 2 dX \y/2J d\ 



which results in 



iB 

\1/ = \I/o = constant and a = «o 1=— X, (79) 

v 2 J 

where ^o, Re(«o) and Im(«o) are constants. These can initially be set to zero in the first cell, 
without loss of generality, but in subsequent cells are generally non-zero. 

Now, the tangent vector in the new cell, k a , is already known from the previous section, and 
i a is the same as in (!76|) , but with n a replaced by n a . The values of \l/o, Re(do) and Im(do) can 
then be found from (!74ll 21 l. Although «o is not needed directly for integrating (!72|) . it is needed 
to find ^o, and so should be recorded. Its value at exit from cell two is found from (!79|) . This 
gives us the required value of ^ — $ when starting the new cell. 



19 They represent the three required degrees of freedom from the eight components of t a , and the five conditions 
it must satisfy. 

20 The symmetry of C a bcd — —C^acd — —C a bdc means that the term involving a in t a does not contribute to 

C. 

21 Equation ([71)1 is 4 equations for 3 variables, so it must be the case that only three of these 4 are independent. 
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The prescription for calculating shear is then the following: 



(1) At the beginning of the cell, take X , fy , <3> and «o f° r that cell. These will have been 
calculated at the end of the previous cell (see steps below). If this is the first cell then take 
X = = % = «o = 0. 

(2) Substitute X , and $o into (1721) . and integrate from the beginning of the cell until 
the trajectory hits another edge. 

(3) Use this X = X(X) to integrate (I6"5j) along the same trajectory, to find r^. 

(4) We now want to calculate Xq and <f>o for the next cell (hats indicate quantities in the 
next cell). Xq is given by Xq = X e , where subscript e means evaluated at the exit from the 
cell. The value of <£>o is given by 

$o = $ e = - arcsin ( ^sin(^ - $o) ) • 



(5) It remains to find and do- These will be found from Eq. ( 1741) . Now, t a (at the exit 
from the first cell) is known from ( 175]) . ( f76l) and ( 1791) . In evaluating t a take B, J, ty , a and 
k a all from the first cell, and use r and A as appropriate at the cell boundary. This gives the 
RHS of flTJ). The LHS is given by 

The t a and k a in this expression are different to i a and k a , and are not just a coordinate 
transformation of them. The tangent vector k is the tangent vector in the new cell (evaluated 
at the boundary still). This is the new vector found in the previous section. The vector i a has 
the same functional form as (176 p . but now with all hatted quantities, as found in the new cell 
(including J from the new cell). 



5.2 Solving the equations. 

The influence of shear complicates the solving of the optical equations (1681) and (1721) . We can, 
however, make quick progress if we are prepared to make approximations. The system we are 
describing involves the propagation of geodesies through multiple lattice cells. At low redshifts, 
and for trajectories that do not pass close by to a central mass, we expect the effect of the shear 
to be small, as the driving term C ~ 1/r 5 . In the absence of shear, when a ~ 0, we can write 
the solution for ta in (165|) as 

r A ~ ci + c 2 A. (80) 

At larger redshifts, and for trajectories that pass close to a central mass, however, the influence 
of shear will accumulate and become non- negligible. 

In a classic paper, Press and Gunn [TTJ addressed the impact of condensed objects on the 
optical properties of the Universe. They modelled the evolution of the optical equations in 
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Figure 10: The top panel shows the evolution of ta as a function of z for a case where o is 
negligible (top, solid line) and for the case where it diverges (bottom, dotted line). In the bottom 
panel we plot the evolution with redshift of the corresponding \a\. 

a clumpy universe as a Markov process, with a performing a random walk through sporadic 
scattering events. It was found that, apart from some extreme scattering events, the cumulative 
effect of the source terms in the Sachs equations is negligible unless the correlation length of 
the objects is very large, or the distances travelled were very great. In fact, the change in shear 
over a Hubble length Hq should be a ~ H and hence a /Ho <C 1. We find that the main 
features derived in [71] are, essentially, in agreement with the numerical results that come from 
integrating the optical equations through a lattice universe. 

In Figure [TD] we plot r^(z) for two different cases; one in which a is always negligible, and one 
in which a increases substantially. Of course, shear is always negligible initially. Far from the 
origin the situation is more complex, and \a\ (or X) can deviate considerably from zero. We find 
that \a\ evolves as if it were subject to a Markov process, performing a random walk and slowly 
diverging away from zero. Once it reaches a critical value, the angular diameter distance, Ta, is 
found to stop growing, and to subsequently decrease, leading to the formation of a caustic. The 
likelihood of this happening for a random trajectory is small, but may become non-negligible at 
moderate to high redshift, when many cells have been traversed. The result, if such events do 
occur, is a substantial suppression of ta, corresponding to a significant brightening of distant 
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sources viewed along the trajectory in question. 

We will now refine Press and Gunn's calculation, and work out how likely it will be for the 
shear to diverge and lead to a turnover in ta- To do so, we first note that there is a high degree 
of regularity and recurrence and that, for the most part, it should be acceptable to consider 
propagation through each domain as almost independent from those that have come before. In 
each domain we rotate the coordinate system so that the photon enters each cell in an upward 
direction (i.e. through the (x,y) plane). By approximating the lattice side as a disc we can 
then write the point of entry into cell % as 

Ti = Lq { = L(qi cos (pi, qi sin (pi, -1/2), 

where L is the side length of the cell, and the probability distributions for fa and (pi are 
P(fa, (pi)dcjid(pi = fadqid(pi, where fa G [0,1/2] and (pi G [0, 2ir]. The direction vector at the 
point of entry can then be written 

Yi ~ fj(sinxcos^,sinxsin^,cosx), 

where x £ [0 ) 7r /2] and £ G [0, 2n]. The vector components in the expressions above are given 
in a Cartesian coordinate system. 

We now wish to find the approximate evolution of X. The X and $ equations tell us that 
$ ^constant, as X becomes large. We also know that \I/ is constant in each domain, but 
rotates randomly between domains. We then have 

N 



X(A)~^Acos(^-^) 



i=0 

where 

Ai = 3mJf / d\^. 

In this limit the trajectory X = X(\) can be seen to become stochastic. 

Ultimately, we are interested in the case where a photon crosses a large number of cells. 
The Central Limit Theorem then implies that X should have a Gaussian distribution. What 
we want to know now is the mean and the variance of this distribution. To find these consider 
that (cos($i - *j)) = and (cos($; - vf^) cos($ j - ^)) = We then have 



(X(A)> = 



1 N 



To proceed further we need to work out an approximate expression for A4. If the angular 
momentum in cell % is given by Jf = (r\ x r{) 2 then we have 

J^ = r^-(r l .r l ) 2 = L 2 ff^-(q,.n,) 2 ]. (81) 
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We now assume that we can approximate each geodesic as a straight line. This gives 



r = Yi + Ti(Xi — X) —n — X'ii 
r 2 = r 2 + f 2 X' 2 - 2r • r.A' 
r 2 = L\ 2 + ffA' 2 -2f,L(q,- ni )A'. 

We will assume that Tj is approximately constant in each domain, and set ta = tai — otiX'. We 
then have that the integral Ai becomes 



A = 3m.] 2 



AAi 



r\dX' 



'o [L 2 q 2 + ffX' 2 - 2f l L(q • n)A'] 5/2 ' 
If we change variable to A = — c^A' then this becomes 



A; 



3mJ 2 



ft.; 



Jo 



(r Ai + A) 2 dA 



LV+(^) 2 A 2 + 2(^)L(q.n)A 



5/2 " 



The further change of variables to Y — ^A then gives 



A, = -3mJ 2 % — 



AYi 



(Y Al + YfdY 



5/2' 



^ 2 7o [(q 2 + Y 2 + 2(q ■ n)F] 

where we have assumed L ^constant over one domain, and defined Ya% = fiTAi/otiL. The upper 
limit of integration is AYi = ~zt' ai ^ i = Ti ^ Xi , which we can approximate to I. The integral 
is therefore time independent, and depends only on the stochastic variables q and n (except 
for the presence of Ya%)- Note, however, that in the limit of multiple scattering events we have 
that Y <C Yai, so that 



AK 



(Y M + YfdY 



'0 [(q 2 + Y 2 + 2(q-n)y] 
We then have that 



5/2 



Y. 



Ai 



AY,, 



dY 



[(g 2 + F 2 + 2(q-n)F] 



5/2 ' 



A = 



3mfia 2 X 2 

-.F(q, n), 



L 2 



(82) 



where 



.F(q, n) 

b 2 - (q • n) 2 ] f 

Jo 



dY 



[(g 2 + F 2 + 2(q-n)F] 5/2 



(q 2 )2 (1 + q ■ n)(2 + 3q 2 + 4q • n - (q • n) 2 ) + q • n(l + q 2 + 2q • n)3((q • n) 2 - 3q 2 



and we have taken ta 
variance of X as 



3(g 2 )i(l + g 2 + 2q-n)% 2 - (q • n) 2 ) 
aX. If we now take fj = fo/a and L = L a then we can write the 

,2-2^4 N 



2L% 



(83) 



i=0 
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The continuum limit of this is 



(X 2 ) 



9m 2 i 2 a 4 (F 2 ) 



pro 



2Ll 



,dr' 



34) 



where we have used Ai = — Ar^/L, and taken the infinitesimal limit so that di = —dr/L. Now 
with (7) = 7/10 we have that a = (r/r ) 2//3 and (3\ = ((r/r ) 22/15 - 1), which allows us to 
perform the integral above to find 



ro \ 4 

Kdr' 



T Q A 

We can then finally write 
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2 9m 2 T 2 a 4 (F 2 }T f, 
\ X > = ccrboi g ( T ) 



22L b p 
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Let us now determine an estimate for the critical value of X at which the shear becomes 
important, and the divergence from the background evolution occurs. First we define the area 
of a bundle of geodesies focused at the observer to be Y = r\. This gives the equation for f^, 

as 

YY _ Y 2 
~Y ~ T 



Now let us take as our benchmark for the point at which divergence begins to be Y = 0. Recall 
that Y ~constant in the absence of shear, and goes from being initially positive (when shear 
is negligible) to being negative when shear causes the turn over and eventual divergence of r&- 
Taking Y = seems like as reasonable place as any to mark the separation of these two regimes. 
At this point we then find the critical value of X to be 



2 



r A r A 




(86) 



where we have used ta = «A, as above (although this is only strictly true when shear is 
completely negligible). The ratio of (X 2 ) to X 2 is then given by 



(X 2 ) 9m 2 i 2 (F 2 )T 



XI 



22Llf3 2 



(t/Tq)" - 1 



2 ' 



*7) 



where (J r2 ) — 10 and Q(t) is given by the expression above. When this ratio becomes greater 
than 1, about 1/3 of the trajectories will begin to diverge^. 



22 Under the assumption of the Central Limit Theorem. 
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Figure 11: T/ie critical redshift, z cr n, as a function of the ratio of mass separation, L , to Hubble 
scale, H . At z = z crit we expect ~ 1/3 of the trajectories will have experienced enough shear 
to cause a caustic. For Milky Way size galaxies L H ~ 10 -4 . 

In Figure [TT1 we plot the critical redshift as a function of the ratio of mass separation, Lq, to 
Hubble scale, Hq. We find that, in the regime that we are considering, only a tiny fraction of 
trajectories should experience enough shear to cause the optical scalars to diverge. For Milky 
Way sized galaxies we find that the critical redshift is z cr u — 2.8, and up until z ~ 1 it should 
be an excellent approximation to take a ~ 0. 

6 Cosmological Observables 

Having discussed the dynamical and optical properties of an Archipelagian Universe in some 
detail, let us now consider the effects such a model of the Universe has on the interpretation of 
cosmological observables. Clearly, the optical properties of the Universe are of great importance 
for very many different types of cosmological observable. These include reconstruction of the 
Hubble diagram, interpretation of anisotropics in the CMB, and relating distances to redshifts 
for the purpose of understanding galaxy number counts and baryon acoustic oscillations. In 
this section we will consider how these types of observables will differ from FRW, in the model 
currently under consideration. 

6.1 Distance measures, and the Hubble diagram 

In this section we will focus on the reconstruction of the Hubble diagram in a universe with 
discrete matter content. In the previous section we argued that shear will play a negligible role 
in the optical equations, at least out to redshifts ~ 1. Very few supernovae have been measured 
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beyond this value, and hence we feel justified in using Equation ( 1801 . 

If we enforce the condition that A = at the observer then c\ = 0, and flSU]) gives oc A. 
We now want to relate A to r. It was found above that 



1 T e T 

1 + Z = — ~ — 
TO \T e 



3 



2<7> +1 M +1 



Integrating this gives 

T e ~^ ~ T ~^ OC A, (89) 

where we have chosen the integration constant such that A — > as r e — > To. This can now be 
substituted into (|45j) to find oc 1 — (1 + z) -1-3 / 2 ^. Etherington's theorem, rr, = (1 + z) 2 rA, 
then gives the luminosity distance as a function of redshift as 

r L = {l + zf - (l + z) 1 -^ (90) 
3 A. _ 2 , (3-2(7)) , , 



- * + V* ~ 4^> J ^(#^ + ° ( } ' (91) 
where the constants of proportionality have been absorbed into Tl- This can be compared to 

r f = z + r Milne = ^ + ^ an d r EdS = ^ + ^ _ + ^4)_ 

For (7) = 1 r£ is somewhere between Milne and EdS, and is actually the same as was 
predicted by Bertotti [SU] . 

r! ertotti = 2: + ii^z 2 + |(g - /)z 3 + {z% (92) 



with the deceleration parameter q = 1/2 and the mass fraction 23 ! / = 0. In this case objects 
at the same z appear dimmer than in FRW. 

For (7) = 2/3 or 7/10 the result is somewhat different. The lower redshift means that 
objects with the same z now appear brighter than their FRW counterparts. More generally, 
from the above we can see that these models have an effective deceleration parameter given by 

C Uce = - 1, (93) 

in the absence of shear. Using our semi-analytically determined value of (7) = 7/10 we find 
that 

r L = {l + zf- (l + z)-? (94) 

~ z -Y/ + f/ + ^- (95) 

As advertised, the luminosity distance as a function of redshift is such that the Universe is 
perceived to be decelerating at an even more efficient rate than in FRW. Indeed, we find 
q = 8/7 for the Archipelagian Universe, as opposed to q = 1/2 for the Einstein-de Sitter 
model. The luminosity distance (|94|) is shown graphically in Figure [T2J in the form of the 
distance modulus^, A dm. 



23 Taking / = 1 in ([§!?)) gives the FRW distances supplied above. 
24 Distancc modulus is denned by Adm = 5 log 10 (r i /r^ Illne ). 
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Figure 12: The distance modulus Adm for the Archipelagian Universe (solid line), an Einstein- 
de Sitter universe (dashed line), a de Sitter universe (dot-dashed line), and for a trajectory with 
(7) = 1 (dotted line). 



6.2 Other observables 



Hubble diagrams are, of course, not the only cosmological observables that are sensitive to the 
optical properties of the Universe. There are various others, and indeed the CMB is famously 
an excellent probe of the dynamics and optical properties of the Universe since the epoch of 
last scattering, before which the Universe was opaque. If we assume the growth of structure in 
the Universe is well specified by the usual perturbative analysis around an FRW background, 
then the principle way in which the CMB will be sensitive to the new distance relations found 
here will be in the relation between angles on the sky today, and length scales at last scattering: 
The projection of the CMB onto our sky. 

A useful probe of such effects on the CMB is the 'shift parameter', which specifies how much 
the peaks in the acoustic spectrum of primary anisotropies are expected to move by. The usual 
convention is to specify this shift with respect to an EdS universe, in which case it is simply 
given by the ratio of angular diameter distances in the test cosmology and the corresponding 
distance in an EdS universe. These distances must be specified so that the Hubble rate at last 
scattering is the same in each of the two universes, so that the physics up until that point is 
also the same. 

From the considerations above we know that the angular diameter distance in the Archipela- 
gian Universe is, in the absence of shear, given by 



lattice 

A ~ 22H l attice 



1 



(96) 
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The corresponding distance in EdS is given by 

r EdS _ 2 



# EdS 



(97) 



As mentioned above, in this kind of procedure one would normally have to relate if EdS to 
jjiattice^ ensur i n g that in both space-times the Hubble rate at last scattering was the same. 
Here, however, we have that the dynamics of the two space-times under consideration are, in 
fact, the same. We therefore have simply that H^ dS = H l Q atttce ^ an d so the shift parameter, S, 
is given by 

SB ;T- 7(i-(i + ,)-T) 



This is shown graphically in Figure 

It is clear that the CMB is very sensitive indeed to the type of modifications to the usual 
distance relations that we are considering. If we are conservative, and only apply (1961) out to 
a distance of z — 2, then we can see from Figure [13] that this already produces a shift in the 
CMB acoustic spectrum of ~ 10%. This effect is already huge, with current experiments able 
to constrain the shift to the order of 1%. 

Extrapolating our results to higher redshifts is somewhat tricky. Of course, we know that 
at distances z > 0(1) we should include the effects of shear. This will likely have a significant 
effect on the angular diameter distance to last scattering. Furthermore, the assumption of Milky 
Way sized central masses for of our cells will likely become a poor approximation at high z, 
especially when we get to the kind of redshifts before these structures even form. As mentioned 
previously, one may consider a different type of lattice in which dark matter particles make up 
the central masses of much smaller lattice cells, although the transition from galaxies to dark 
matter particles will certainly be a highly non-trivial process. 

For these reasons our numerical codes are currently unable to make predictions out to the 
very high redshifts required to perform a proper analysis of the CMB, but the discussion above 
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indicates that there could be significant deviations from the standard FRW picture. We will 
leave further considerations about CMB observations in these cosmologies, which will likely be 
a highly non-trivial matter, to future publications. 

As well as Hubble diagrams and the CMB there are a wealth of other observables available 
to cosmologists. Many of these, such as baryon acoustic oscillations and galaxy number counts 
are highly sensitive to the relation between cosmological time and redshift. In our models, this 
relation can be read off from Equation f)56p . The smaller redshifts, that appear to be a generic 
prediction of this model, mean that baryon acoustic oscillations should be measured at a lower 
redshift in these models, as compared to the case of FRW. In terms of galaxy number counts, 
the lower redshift will mean that more spatial volume will be included out to some particular 
z, and hence one should expect a higher density of galaxies. 

There are clearly many observational consequences of treating the matter distribution as 
being discrete, rather than continuous. It is not the purpose of the current paper to work 
through all of these rigorously, or to provide a complete working alternative to the standard 
ACDM concordance model. Rather, we have aimed at gaining a decent understanding of the 
simplest dynamical and optical properties of the simplest models. Once again, more thorough 
investigations of the generalisations of these models, and their observational consequences, will 
be postponed to future studies. 

7 Discussion 

In this paper we have considered an inhomogeneous cosmological model with basic properties 
that are, in essence, similar to our own Universe. The idea under-lying the model is that we 
should break away from the doctrine that the Universe is permeated by a continuous fluid of 
matter, and instead consider the view that it is predominantly empty space, punctuated by 
islands of mass. We identify these islands as the de facto building blocks of structure: Galaxies 
with a mass similar to that of the Milky Way. 

Our study extends and applies the ideas of the lattice universe model proposed by Lindquist 
and Wheeler [52] that has, until now, remained relatively unexplored. We expect this approach 
to cosmology to allow fresh insights into the effects of inhomogeneity in the Universe, and, 
perhaps, to allow fresh approaches to understanding some of the unresolved problems that 
have recently been revealed by observations. The most notable of these is the evidence for Dark 
Energy. Such evidence relies heavily on understanding how to relate luminosities and redshifts 
of distant astrophysical sources to the expansion history of the Universe, and although we do 
not uncover any evidence that Dark Energy itself could be explained by considerations of the 
type explored in this study, this does not mean that refinements of the usual FRW cosmology 
will not be of use for properly interpreting observational data^l. Neither does it mean that 
future generalisations and refinements of this model will not be more successful in this regard. 
We consider the exploration of the optical properties of different space-times to be essential to 

25 In [HI] we use Schwarzschild-de Sitter cells, instead of Schwarzschild cells, to construct a lattice universe. In 
this way it is be possible to study the effects of discretization of the matter content on estimates of cosmological 
parameters in a more realistic fashion, by including the effects of A. 
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gaining a full understanding of observations in cosmology. 

Although we have chosen to study a highly symmetric model - evenly spaced masses with 
a critical density - we believe that it does allow us to reach some conclusions that should have 
consequences for the real Universe. For a start, the structure of the Sachs optical equations 
are such that a bundle of null rays travelling through empty space, and with negligible shear, 
should expand to give an angular diameter distance ta oc A. This corresponds to sources at the 
same affine distance being dimmer in the Archipelagian Cosmology than they would be in the 
corresponding Einstein-de Sitter universe. It appears to us that this should be a robust generic 
prediction for any space-time which is mostly empty space. If the photons themselves do not 
pass through the notional continuous fluid, then they do not experience the extra focusing that 
such a fluid would produce. Photons experience the integrated effects of the geometry through 
which they pass, and not the 'average' global geometry. 

Similar effects have also been found in the context of weak lensing. In most studies on this 
subject it is usual to take the approach used by Press & Gunn [TTJ whereby one considers a 
tube of FRW space-time, removes the evenly distributed mass, and replaces it with the equiva- 
lent density of compact objects. The cosmological dynamics are then assume to remain FRW, 
while for the purposes of integrating the optical equations one assumes a vacuum. It has been 
shown by Weinberg [72] and Kibble & Lieu [73] that this approach leads to a distribution of 
magnitudes with the same average as would have resulted from the continuous energy den- 
sity^. However, the distribution of magnitudes is highly skewed [75], with a large number of 
trajectories experiencing little shear, and a small number experiencing large shear due to the 
occurrence of rare lensing events. The sample mean of a limited number of supernovae at small 
z can therefore result in deviations from the FRW average [76], and the skew can be used as a 
test of the discreteness of the matter content [TTJ . Although our model is different to those just 
mentioned, our results with regard to the effect of lensing do not appear to be in disagreement. 

We also find, however, that redshifts can be altered from their FRW values. We have shown 
that, in the context of our model, the redshift can be reasonably well understood in terms of the 
relative motion of the boundaries of the lattice cells. Unlike in an FRW cosmology, however, it 
is only radial geodesies that experience redshifts that would correspond to the Doppler shift due 
to the recessional velocity of the global expansion. All non-radial geodesies experience smaller 
redshifts, as they see a locally anisotropic geometry, and not the global 'average'. Despite the 
additional complications that arise from this anisotropy, we find that after traversing a large 
number of cells typical trajectories quickly approach a mean value 1 + z ~ (1 + zfrw) 1 ® ■ 

These results allow us to then calculate observable measures of distance as function of the 
redshift, as would be recorded by astronomers. We find that although objects appear dimmer at 
the same affine distance, they appear brighter at the same redshift. This is due to the redshift 
being lower to distant objects, as discussed above. Hence, we find that the distance moduli 
measured in a critical density Archipelagian Universe should fare even worse when compared 
to supernova observations than a standard Einstein-de Sitter universe. This is so even in the 
case of negligible shear, and can be made worse when shear is included. Other cosmological 
observables are also discussed, and the influence of the current considerations on interpreting 
them is outlined. 

26 See, however, [71) who point out the occurrence of caustics changes this result. 
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An interesting effect that does not occur in FRW cosmology, but seems to be an inevitable 
feature of Archipelagian Cosmology, is the formation of caustics in the trajectories of photons. 
The discrete nature of the mass in this model means that shear, even if it is negligible initially, 
will eventually accumulate, leading to the focusing, and hence extreme brightening, of some 
distant sources. We find this effect to be small out to redshifts of a few if the masses involved 
are sufficiently smalQ but that it increases substantial if we consider large agglomerations 
of mass, and large redshifts. Such considerations may therefore be important when observing 
distant objects, and for CMB observations, in particular (see [71]). 

The Lindquist- Wheeler model is highly idealised, and is only an approximate solution of 
Einstein's equations. Nevertheless, it has allowed us to draw conclusions about space-times 
which are qualitatively different from the standard FRW space-time. An obvious next step to 
improve the realism of the model is to consider irregular lattices, such as those provided by, 
for example, a Voronoi tessellation. This is an ambitious proposal, and will undoubtedly throw 
up many difficulties that the Lindquist- Wheeler model avoids. However, it could allow us to 
model the observed Universe in a new way, and provide a new laboratory to study, for example, 
non-perturbative back-reaction. Improving on the accuracy of the geometry as a solution of 
Einstein's equations could be more difficult still, and the only way we can see to make progress 
on this front is by using either numerically methods or weak field approximations. 

Finally, we should be careful to ensure that systematic errors have not unduly influenced 
our results. While it is routine to assume that the vacuum space-time geometry outside mass 
concentrations such as galaxies is well approximated by the Schwarzschild solution, this is not 
exactly true: There will be external influences from other gravitational sources that will deviate 
from perfect spherical symmetry. Determining the extent to which such influences effect the 
result we have found here is likely to be a complicated matter, and will be the subject of future 
publications. One may also be concerned as to the possibility of introducing systematic errors 
from the matching conditions at cells boundaries. Here we have used two different methods 
to propagate trajectories between cells, and found that our results are largely insensitive to 
which is chosen. This suggests that the boundary conditions are not a source of considerable 
error. Furthermore, we have found analytic approximations that are in good agreement with 
our numerical integrations, and that allow us insight into the physical origin of the effects we 
have found. We also show in a separate publication that the familiar optical relations of de 
Sitter space are recovered in the limit — > 1 [61]. We consider these results to suggest that 
the systematic errors in our model are indeed under control. Confirming this by determining 
the observational consequences of the approximations used in this model will be the subject of 
future work. 

In conclusion, we have shown that the differences between making observations in a standard 
Einstein-de Sitter universe, and in a universe with discretized matter content, can be substantial. 
A lesson we take away from this study is that the assumption of continuity of energy in the 
Universe, while appearing innocuous, can affect our interpretation of redshifts and luminosity 
distances considerably. In order to make precise statements about the Universe, we should 
therefore make sure we take these issues into account. 

27 For example, if they are Milky Way sized. 
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A Tiling the 3-space 

In order to consider a lattice model of the Universe, it is first necessary to consider how to build 
a lattice out of tessellating cells that fill the 3-space. This problem is known as 'tiling', and has 
been considered in detail by Coxeter [64] . 
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Table 1: Polyhedra that tile 3-surfaces of constant curvature, the number of cells required to fill 
the space, and the structure of the lattice (in the form {pqr}). Hyper-spherical space is denoted 
by '+', flat 3-space by '0 7 , and hyperbolic 3-space by '— See the text for an explanation of 
{pqr}. For further details see [64] . 

In Table [TJ below, we show all the regular polyhedra that tile all 3-spaces of constant 
curvature. For the case of a positively curved 3-space, as considered by LW, there are 6 
possible tilings with regular polyhedra. These are with N=5, 8, 16, 24, 120 and 600 cells. For 
flat 3-space there is only one possible tiling, consisting of an infinite number of cubes. Lastly, 
for negatively curved, hyperbolic 3-space, there are 4 polyhedra that can completely tile the 
space. Lattices constructed on flat and negatively curved backgrounds have been considered in 
a cosmological context by Redmount in |78j . 

The structure of a lattice can be well described using the compact notation {pqr}, which 
is given for the possible tilings of 3-spaces of constant curvature in Table [H In this notation p 
denotes the number of edges on the face of a lattice cell, q denotes the number of faces that 
meet at the apex of any individual cell, and r denotes the number of cells that meet at an 
edge. Hence we have {434} for the cubic tiling of a flat space, as there are p = 4 edges to the 
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square face of a cubic cell, r = 3 squares meeting at the corner of each individual cell (if it 
were considered in isolation from the other cells), and r = 4 cubes meeting around every edge 
of every square. 

For discussion of the efficacy of replacing the above polyhedra with spheres, the reader is 
referred to |55|. 



B Geometric Set-up With a Cosmological Constant 

It is straightforward to generalise the lattice cosmology to include a non-zero A. The metric in 
each lattice cell simply becomes Schwarzschild-de Sitter, with the line-element 

ds 2 = -(l- — - -A dt 2 + 7 . r + rW. (99) 

A coordinate system in which space-like surfaces, with r =constant, are orthogonal to the 
world-lines of elements of radially free-falling time-like shells is then given by the coordinate 
transformation (with E > 1) 

dr (E — (l — — — ±r 2 )) dr , 
dt = | i i r - ^ — (100) 

Ve (i-^_A r2 ) y j E 2_ E ^ 1 _^_A r2 y 

so that the line-element fl99|) becomes 



, 2 If 2m A 2 \ , 2 2 / / 2m A A , , dr 2 2 ~ , . 
ds 2 = -—(l r 2 dr 2 - —JE - 1 r 2 jdrdr + — + r 2 dVl 2 , (101) 



E \ r 3 J E y V r 3 / E 

Under the re-labelling r — > a, m — > M/2 and E — > 1 — k we then have the equation of motion 
for a time-like particle in radial free-fall being given by 

a 2 M k A , . 

+ p (102) 

which is the Friedmann equation for dust, with a cosmological constant. 
In this case the normalised four-velocity of a free-falling object is 



u a = \l H lE-(l- ^-^,0,0j, (103) 

and it can be seen that such trajectories are orthogonal to surfaces of constant r as 

u a n a = (104) 
for any arbitrary vector n a = (0; n r , n 9 , n^) that exists in such a surface. 
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C Geodesic Equations in Cartesian Coordinates 



For E — 1 we can also express the geodesic equations in terms of Cartesian coordinates, rather 
than spherical polars. The line-element ([8]) then appears as 



(2m\ 2m 
1 Jdr 2 - 2d—(xdx + ydy + zdz)dr + dx 2 + dy 2 + dz 2 . (105) 

In terms of these variables the Euler-Lagrange equations read 



2m\ . 2m t 



l- — )i+\/—(xx + yy + zz) = B (106) 



x y z 2m.. m. 2 

- = - = - = \ —t-—t\ (107) 

x y z V ?" r d 

with the null constraint 

^(xi + yy + zi) 2 + ( 1 - — ] (x 2 + y 2 + z 2 ) = B 2 . (108) 

r i V r / 

These equations do not allow integrals as easily as in the spherical polar case, but have fewer 
problems with coordinate singularities. 



D Geodesic Equations With a Cosmological Constant 

The Euler-Lagrange equations for null geodesies in the space-time specified by (11011) can be 
written 



d f 2m A 9 \ / / 2m A , . , 

1 r 2 ]t + a E- 1 r 2 r = (109) 



dX 



dr / / 2m A A df ( m A \ , ^ ao „ 9 ^ :o , 

— -\E- 1 r 2 ) — + I — r f 2 = £r# 2 + £r sin 2 6<p 2 110 

dX \ \ r 3 J dX \r 2 3 / 

(r 2 0} = r 2 sm9cos9<p 2 (111) 

Vsin 2 #0) = 0, (112) 

HA 

together with the null constraint 



d 

dX 

d 

dX 



1 - ^ - f 2 + f 2 + Er 2 6 2 + Er 2 sin 2 6<p 2 - 2 J E - (l - ^ - ^r 2 ) rf = 0. (113) 
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Figure 14: The spherical boundary of a cell. We can consider the photon to enter at the South 
pole, without loss of generality. 



E Finding (7) 

To find (7) consider a photon entering the spherical cell, rather than the cube, as illustrated in 
Figures [TH and [J31 

Now consider a random set of trajectories in flat 3 dimensional space (random refers to the 
measure of trajectories which is invariant under translations, rotations and reflections). These 
trajectories represent the paths of photons in our model. We now want to know the distribution 
of chords that result from these trajectories intersecting a unit sphere that is placed in this space. 
It is shown in [79] that this situation is exactly equivalent to selecting at random two points on 
the surface of the sphere and joining them with a chord. In this case we can choose a coordinate 
system so that the first point is at the South pole. For the chord in question the ratio a./j3 is 
then given by xjl — cos 6, from Figure [TBI We will now proceed as prescribed by Berengut [79] 
to find the mean and variance of cos#. 

The probability of finding our second point a chord length I < L < I + dl away is given by 
the area of the band in Figure [TH over 4n (we consider a unit sphere for now) . To find the area, 
A, of this band consider the cartoon schematic shown in Figure [T5j It can be seen that 



dl = dS cos = dS sin 9, 



(114) 



as <p = 7r/2 — 9. The band area is then 



A = dS x 2nl sin 9 = 2irldl = 8ir cos 9d cos 9 



(115) 



as / = 2 cos#, and the probability of finding our second point in this band is 
P(cos 9 < cos < cos 9 + d cos 9) = — = 2 cos 9d cos 9, 



(116) 
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giving the distribution /(cos#) = 2cos6>. The mean and variance of cos 9 are then given by 
cos 9 = 2/3 and a^ os9 = 1/18. 

The value of (7) will then be given by the mean of the sample of n cells that the trajectory 
passes through, giving 

(7) = I (117) 

with variance 

For a large number of cells this variance will soon become negligible, and we will have (7) — » 2/3, 
as n — > 00. 
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